Systematic Derivation of Order Parameters through Reduced Density Matrices 
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A systematic method for determining order parameters for quantum many-body systems on lat- 
tices is developed by utilizing reduced density matrices. This method allows one to extract the 
order parameter directly from the wave functions of the degenerate ground states without aid of 
empirical knowledge, and thus opens a way to explore unknown exotic orders. The applicability of 
this method is demonstrated numerically or rigorously in models which are considered to exhibit 
dimer, scalar chiral, and topological orders. 
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Determining order parameters is one of the most im- 
portant issues in the study of many-body systems. A 
suitably chosen order parameter for a symmetry-breaking 
phase provides an intuitive picture of the long range or- 
der and is the necessary starting point of the Landau- 
Ginzburg-type effective description of the system 0. 
Combined with Wilson's idea of renormalization group 
Q, such an effective theory becomes a powerful tool in 
analyzing the nature of the phase transitions to other 
phases. 

In spite of the importance of this issue, a general 
method to obtain an order parameter in a given model is 
not available. The knowledge of previous examples sug- 
gests some candidates but this empirical method may fail 
in the case of a new order. Especially in a system with 
strong frustration and/or quantum fluctuation, the order 
parameter can be quite non-trivial. With more examples 
of exotic orders becoming a subject of great theoretical 
and experimental interest, a systematic method for de- 
termining an order parameter would be strongly desired. 

In this letter, we present a solution to the quantum ver- 
sion of this problem. In the quantum case, when a sym- 
metry of the Hamiltonian is broken spontaneously in the 
thermodynamic limit, there appear degenerate ground 
states (GS). An order parameter can be identified with 
an operator which distinguishes the degenerate GSs. The 
central idea of our method is to search such an operator 
by comparing the reduced density matrices (RDM) of 
the degenerate GSs for various subareas of the system. 
A RDM efficiently encapsulates the expectation values 
of all the operators on the concerned area. If the RDMs 
of the GSs are different on an area, an order parameter 
can be defined on that area. In this way, we can de- 
termine the smallest area on which an order parameter 
can be defined. Moreover, for the resultant area, we can 
construct an "optimal" order parameter from the RDMs. 
This method can be applied to the low-energy eigenstates 
obtained by exact diagonalization, for instance, and can 
reveal the order parameter without bias. 

We will demonstrate the effectiveness of this approach 
in concrete models. We will consider the multiple-spin 



exchange model on the ladder and detect dimer and 
scalar chiral orders which have been found in previous 
studies 0, . We will also consider a resonating valence 
bond (RVB) liquid in a solvable quantum dimer model 
(QDM) 5] and rigorously show that its GSs cannot be 
characterized by any local order parameter |l'2^ . Namely, 
the model will unambiguously be shown to posses a topo- 
logical order. 

Methodology — Suppose that we have obtained the 
low-energy spectrum and eigenstates of finite-size sys- 
tems by exact diagonalization, for instance. In a phase 
breaking a discrete symmetry, we find a finite number of 
nearly-degenerate GSs which become asymptotically de- 
generate when increasing the system size. Each of these 
states does not break any symmetry of the Hamiltonian 
but its quantum numbers indicate what symmetries are 
broken in the thermodynamic limit. 

Let us focus on the simplest case: the Hamiltonian is 
invariant under the translation by one lattice spacing, 
T, is real in terms of {SJI-basis, and exhibits doubly- 
degenerate GSs, |$i) and |*2), with momenta k = and 
7T, respectively. In this case, we expect the breaking of 
the translational symmetry (doubling of the unit cell) in 
the thermodynamic limit. We set |*i) and |*2) real, 
i.e., JC\$i) = (i = 1,2), where /C denotes the time- 
reversal operator which converts every component of a 
vector into its complex conjugate in terms of {S'|}-basis 
0. 

We construct the symmetry-breaking GSs, |*i) and 
|*2), a s linear combinations of |*i) and |*2)- We re- 
quire that they be orthogonal ((*i|*2) = 0) and be 
exchanged under T (T|*i( 2 )) oc |*2(i)))- There are 
still two possibilities, depending on whether the time- 
reversal symmetry is broken (/C|*i( 2 )) oc |*2(i))) or n °t 
(/CI* 1(2)) oc 1*1(2)))- For each case, the symmetry- 
breaking GSs are constructed as 

1*1,2) = (|*i) ± |* 2 ))/v / 2 (/C-unbreaking case), 
|*i,2) = (|$i)±i|# a »/V2 (^-breaking case). (1) 

Here the both possibilities have to be examined since, 
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due to the anti-unitarity of the time-reversal operator, 
we do not know from the quantum numbers whether the 
system breaks the time-reversal symmetry or not. 

Next we search an operator which distinguishes the 
symmetry-breaking GSs by comparing the RDMs p l n — 
Ttq (i = 1, 2), where f2 is an area in the system 

and Cl its complement. To quantify to what extent the 
two RDMs are different, we introduce a measure as 



max 

|O n |<l 



Tr(OaP^) 



(2) 



where On is a variational (hermitian) operator on Q sat- 
isfying I {ip I Oq, < 1 for any normalized vector \ip). If 
diff /Jq) > 0, there exists an operator on Q, distin- 
guishing and This measure has the following 
useful properties, a) normalization to a definite range 
< diff(j0f2, Pq) < 2, for an arbitrary area fi. b) mono- 
tonicity: if an area A completely contains an area fl, we 
have diS(p Q) p Q ) < dffi(p\,f%). 

Using the eigenvalues {Xj} and the eigenvectors {\j}} 
of Apn = p^ — Pq, Eq. (J2J) can be simplified as 
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FIG. 1: Numbering of the sites on the two-leg ladder. 



Area Q 


9 = 0.12tt 


9 = 0.19tt 


diffl diff2 


diffl difT2 


{1} 
{1,2} 

{1,1'} 
{1,2'} 
{1,2,1'} 
{1,2,3} 
{1,2,1', 2'} 


0* 0* 
0.5698 0* 
0* 0* 
0* 0* 
0.5698 0.0267 
0.6579 0.0670 
0.5698 0.0462 


0* 0* 
0.0029 0* 
0* 0* 
0* 0* 
0.0029 0.3340 
0.0033 0.2365 
0.0029 0.5785 



TABLE I: Values of "diff" for various areas in the 14 x 2 
ladder. The points, 9 = 0.12-7T and 0.197T, are the representa- 
tive points respectively in the staggered dimer and the scalar 
chiral phases found previously The sites are numbered 

as shown in Fig. Some zeros (indicated by *) are exact 
consequences of the symmetries. 



diS( p n ,p Q ) = max 
|O n |<l 
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^{j\On\j) 



(3) 



Here the maximization is done by the "optimal" order 
parameter: 



O, 



(opt) 



= E \i) s § n a j 01' 



(4) 



where sgn Xj is the sign of Xj if Xj 7^ and is zero if 
Xj = 0. Both the measure and the optimal order parame- 
ter can be calculated by (numerically) diagonalizing Ap^. 
As we have discussed above, generally we have to exam- 
ine both the /C-unbreaking and /C-breaking combinations 
in eq. Q. In the following, we denote the measure (0 
for the /C-unbreaking and breaking cases as "diffl" and 
"diff2" , respectively. 

The generalization of this method to systems with 
more than two degenerate GSs can be formulated as an 
optimization problem, which will be presented elsewhere. 

Simple examples — To illustrate this method, let us 
consider two simple examples, Neel and dimer orders. 
The corresponding symmetry-unbreaking GSs are respec- 
tively given by 
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U---)±UT •••}), 

\sn-i,n) ± 1^23} 



(5) 

|sjv,i)),(6) 



where |sy) denotes a singlet bond. For the Neel order, 
diffl = 2 for a 1-site area {1} while diff2 is zero for the 
same area. Thus the optimal order parameter should 
be constructed from the RDMs of the /C-unbreaking GSs 



on {1}, resulting in C'|°j t ' ) = 2Sj\ For the dimer order, 
on the other hand, both "diffl" and "diff2" are zero for 
{1} but we find diffl = 3/2 > for a 2-site area {1,2}. 
The resultant optimal order parameter Q for this area is 



O 



(opt) 
{1,2} 



2Si ■ S2 + 1/2. We have obtained the expected 



order parameters for both of the simple examples. 

Actually, in order to establish the presence (or absence) 
of an order parameter on a given finite area, the mea- 
sure "diff" @ has to be defined in the thermodynamic 
limit. However, in most of the applications, especially 
in numerical calculations, we would be only able to cal- 
culate the corresponding quantity in finite systems. We 
expect that, in a gapped system, the "diff" should con- 
verge exponentially to the true value, when the system 
size is taken to infinity. While the systematic study of 
such finite-size effect is outside the scope of the present 
Letter, the following application demonstrates the use- 
fulness of our measure even in numerical diagonalization 
of relatively small systems. 

Application I — We consider the 2-spin and 4-spin ex- 
change model with spin S = 1/2 on the two- leg ladder: 



H 



(7) 



where the two summations run over the (vertical and hor- 
izontal) bonds and the squares, respectively. According 
to recent analyses 0, Q , two ordered phases breaking the 
translational symmetry have been found: the staggered 
dimer phase (0.077T < 9 < 6 C ) and the scalar chiral phase 
(&c < S < 0.397t), separated by the exact self-dual point 
9 C = tan- 1 (l/2) ~ 0.1476tt 0; see Fig. ^a). 
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FIG. 2: (a) Phase diagram obtained from earlier studies 
HQ. (b) Values of "diffl" and "diff2" for the fixed ar- 
eas, {1, 2, 1'} and {1, 2, 1', 2'}, versus 0. The values of "diffl" 
(filled symbols) for the two areas are exactly the same. As for 
"diff2" (open symbols), the upper and the lower points refer 
to {1, 2, 1', 2'} and {1, 2, 1'}, respectively. Lines are guides to 
eyes. The vertical dotted straight line represents the self-dual 
point 9 C . Our calculation assumes a quasi-degeneracy of the 
GSs and is therefore invalid for 9 < 0.07-7T and 0.39-7T < 9. 



In both regions, the finite-size spectra obtained from 
exact diagonalization exhibit two nearly-degenerate sin- 
glet GSs with quantum numbers, (k Xl k y .a) — (0,0,1) 
and (7r,7r,— 1), where (k x ,k y ) denotes the momentum 
and a the reflection with respect to a rung. We con- 
structed symmetry-breaking GSs from these states and 
calculated "diff" for various areas; see Table Q] At 
9 = 0.12-7T, the minimum area required to find an or- 
der parameter is {1,2} and the time-reversal symme- 
try is unbroken. Since Ap^i ^} is symmetric under the 
spin rotations and the time reversal, it must be propor- 
tional to S± ■ S 2 and hence the optimal order parameter 
is 0j° p *j = 2Si ■ S 2 + 1/2. At 9 = 0.19vr, the minimum 
area consists of three sites (e.g., {1, 2, 1'}) and the time- 
reversal symmetry is broken. Since Ap{ 12 ,i'} is symmet- 
ric under the spin rotations and antisymmetric under the 
time reversal, we have Ap^ 2 ,v] oc S\ ■ (S 2 x Si>), and 

hence O^i^i'} — -^Si ■ (S 2 x Si'). In this way, we have 
derived the dimer and the scalar chiral operators as the 
appropriate order parameters in a systematic way. 

In Fig. Ufa), (^-dependence of "diffl" and "diff2" are 
shown for fixed areas. Rapid changes in the values of 
"diff" can be seen around the self-dual point 9 C , confirm- 
ing the phase transition between the two ordered phases. 
For {1,2,1', 2'}, the values of "diffl" and "difT2" cross 
exactly at 9 C . For {1,2,1'}, the crossing of "diffl" and 
"diff2" deviates from 9 C but approaches it when increas- 
ing the system size. In general, such a crossing indicates 
a phase transition between ordered phases which cannot 
be distinguished by the GS quantum numbers. 



Application II — We next consider a solvable QDM on 
the kagome lattice introduced recently 5] (for a review, 
see also section 5 of Ref. 0), which is one of the mi- 
croscopic models realizing a short-ranged (so-called Z 2 ) 
RVB liquid. This model also provides an example of a 
solvable Hamiltonian |jj for a topological quantum-bit 
based on a QDM Before applying our method, we 
briefly review the definition of this model and some basic 
concepts. 

This model is simply expressed in terms of the arrow 
representation [|J of dimer coverings defined in the fol- 
lowing way. The sites of the kagome lattice K can be 
identified with the centers of the bonds of the hexagonal 
lattice H. For a dimer covering of K, we assign orien- 
tations (arrows) to the bonds of H so that the arrow on 
each site of K points towards the interior of the trian- 
gle of K where the dimer occupying the site is (see Fig. 
Ota)). As a consequence, the number of incoming arrows 
is even (0 or 2) at every triangle. Let S be the set of ar- 
row configurations satisfying this local parity constraint 
at every triangle. There is a one-to-one correspondence 
between S and the set of all dimer coverings. 

We define r z (i) as the operator which flips the arrow 
on the site i of K. Dimer movements can be represented 
as loop products of r 2 operators. The Hamiltonian we 
consider is the sum of the loop operators around the 
hexagons h of K: H = - J2 h 111=1 ^(ih^), where i h>a 
are the six sites of the hexagon h. 

If this model is defined on a surface with a non-trivial 
topology (cylinder, torus, etc.), arrow configurations in 
S can be grouped into topological sectors which are not 
mixed by any succession of local dimer moves. From 
now on, we concentrate on the case of the cylinder for 
simplicity (but all the results can be easily generalized 
to other topologies). We draw a cut A (passing through 
the bonds of H) going from the top to the bottom of 
the cylinder. We classify arrow configurations into two 
topological sectors, S + and S~ , depending on whether 
the number of arrows crossing A to the right is even or 
odd. The spectrum can be determined separately in each 
sector. Using the standard Rokhsar-Kivelson argument 
|lfj| , one can show that the ground state in a given sector 
is exactly the equal-amplitude superposition of all dimer 
coverings (arrow configurations) belonging to that sector: 



I/O 



E i«>. 



(8) 



aESi 1 



These two states are exactly degenerate and form a 2- 
dimensional GS subspace. 

Now we consider the RDM of a state in the GS 
subspace and discuss how it depends on the choice of 
The area ft is given as a set of bonds of H. The RDM 
pa = Trj2 |^ r )(^| is defined by tracing out the degrees of 
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FIG. 3: (a) Dimer covering and arrow representation. 
(b)Local area Q, on the cylinder. The cut A and the loop 
A* can be drawn so as not to touch Q,. 

freedom (arrows) on 57: 

(ai|pn|oa) = 5^<oi,o|*)(*|o2,o), (9) 

a 

where a\ and ai are arrow configurations on Q and the 
sum is over all the arrow configurations a on fi. By ex- 
pressing |*) = X^^mIm) with \a^\ 2 = 1, pu can be 
expanded as 

pn = Y, a »< M n, M% , = Tr\Mv\. (10) 

First we assume that is a (finite) local area ; see Fig. 
Gib). We prove the following relations: 

M+~=0, M% = M%, (if) 

where the overbars represent the sign flip. To prove the 
first relation, we choose the cut A so as not to touch f2. 
Then the parity along A for an arrow configuration de- 
pends only on its part on f2. Let us consider the matrix 
element of Mq": (ai\M^\a 2 ) = S g {ai, a\+)(— |a2, a). 
Since the two configurations, (ai,a) and (a2,a), have 
common parity, (ai,a\+) and (—\a2,a) cannot be non- 
zero at the same time and hence we obtain Mq^ = 0. 

To prove the second relation in Eq. i|ll[l . we draw a 
loop A* (passing along the bonds of H) encircling the 
cylinder so as not to touch the area fL We introduce a 
loop operator along A*: T£„ = IL s a*i~ z (?). This opera- 
tor acts only on f2 and maps |±) to showing Mq" = 
Tro (TZ.\p)(v\TZ.) = Ttq (\p)(D\(T^) = Mg. 

From Eqs. (|f 0|l and (fill) , we see that pn is indepen- 
dent of the choice of \^f) . Thus no order parameter can be 
defined on an arbitrary local area. In the context of topo- 
logical quantum bit based on QDMs |3,|8j, this shows the 
stability of quantum information against external noises 
(coupling locally to dimers). 

The situation is different if f2 has a non-trivial topol- 
ogy, namely, extends from the top to the bottom of the 
cylinder (case A)) or encircles the cylinder (case B)). In 
case A), we can choose the cut A inside f2. Then the 
parity along A can be expressed as the operator acting 
on £1 and distinguishes the different topological sectors. 
It can be considered as a non-local order parameter dis- 
tinguishing |+) and |— ). Similarly, in case B), two states 
(|+) ± |— ))/v2 are distinguished by the loop operator 
T?» with A* defined inside f2. 



We have shown that the GSs of this QDM cannot be 
distinguished by any local operator but by non-local op- 
erators defined on areas with non-trivial topologies. We 
comment that a similar result has been shown without 
using RDMs by Ioffe and Feigel'man in their study 
on a related model. We stress, however, that our formu- 
lation based on RDMs has an advantage in its generality. 
As demonstrated in the ladder model, it can be applied 
to models without exact solutions, by combining it with 
numerical calculation, for example. 

Conclusions — We have developed a method which can 
determine order parameters without using any empirical 
knowledge. The two applications confirmed its applica- 
bility to exotic orders and, especially, its relevance for an- 
alyzing topological orders. We expect that our method 
will shed some light on the controversies in some frus- 
trated quantum antiferromagnets (see Ref. [(| and ref- 
erences therein), e.g., the J\ — J2 model on the square 
lattice and the multiple-spin exchange model on the tri- 
angular lattice. 
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